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PROGRESS  REPORT  ON 

'MODELING  AND  STATISTICAL  ANALYSIS  OF  BIOASSAY  DATA" 


Donald  P.  Gaver 
Professor  of  Operations  Research 
Naval  Postgraduate  School 
Monterey,  CA  93943 


1.  BACKGROUND 

The  objectives  of  the  above  project  were  formulated  in  discussion  with  Mr. 
Henry  Gardner  of  U.S.  Army  Medical  R&D  Command,  Ft.  Detrick,  MD.  The 
project  purpose  and  workscope  was  stated  in  a  MEPR,  effective  Aug.  05, 1991,  as 
follows; 

The  Department  of  Operations  Research,  Naval  Postgraduate  School, 
Monterey,  California,  shall  conduct  quantitative  analysis  of 
environmental  bioassay  data  to  be  provided  by  the  U.S.  Army  Biomedical 
Research  and  Development  Laboratory  (USABRDL).  Individual  cancer 
bioassays  shall  be  evaluated  as  well  as  all  field  test  data.  Meta-analytical 
methods  for  assessing  integrated  biological  assessment  data  shall  also  be 
accomplished. 


2.  APPROACHES  TAKEN,  AND  PROGRESS 
a)  Data  Availability 

Detailed  data,  concerning  pathologist's  assessments  of  physical  changes  in 
(medaka)  fbh  organs  possibly  resulting  from  toxin  dosage,  were  sent  to  the 
investigator  on  Dec.  19,  1991;  they  arrived  about  a  week  later.  Summarized 
data  were  available  somewhat  earlier,  i.e.  at  the  end  of  November.  The  form  of 
the  data  clearly  influences  the  type  of  statistical  analysis  and  modeling,  so  our 
effort  was  directed  towards  formulating  appropriate  model  types  and 
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developing  techniques  for  applying  or  fitting  these  to  data.  Actual  applications, 
i.e.  fits  amd  error  assessments,  plus  comments  about  the  apparent  sensitivity  of 
results  to  the  experimental /test  designs  (numbers  of  subjects /fish  in  groups, 
intervals  between  sacrifice,  etc.)  have  been  made  in  Appendix  C  but  require 
extension.  They  can  be  enhanced  during  a  subsequent  contract  phase,  which 
will  be  proposed. 

Id  Model  Development 

The  objective  of  the  project,  to  date,  has  been  to  focus  on  biological  issues 
believed  to  be  important  in  converting  toxin  dosage  to  pre-cancerous  and 
cancerous  cells,  and  to  translate  these  into  appropriate  quantitative 
mathematical  terms.  In  particular,  models  that  stem  from,  the  clonal 
expansion  mechanisms  identified  by  Moolgavkar  and  co-workers,  cf.  sample 
references  (1979),  (1983),  have  been  studied,  generalized  (to  account  for  possible 
variability  or  susceptibility  between  individual  fish),  and  adapted  for  fitting  to 
actual  data.  For  an  example  see  Appendix  C,  where  numerical  results, 
including  uncertainty  analyses,  are  presented, 
c)  Some  Details  Concerning  Model  Conception 

It  appears  to  be  widely  believed  that  pre-cancerous  conditions  in  an  organ 
(the  liver)  occur  as  a  result  of  cell  clonal  expansion,  followed  by  a  promotion 
(to  tumor)  event.  Specific  models  for  this  has  been  proposed  and  developed  by 
Moolgavkar  and  co-workers.  More  recent  work  is  by  C.  }.  Fortier  and  co¬ 
workers.  References  appear  later. 

The  basic  mechanism  is  treated  as  random  or  probabilistic:  an  initiating 
event,  e.g.,  caused  by  contact  with  toxin,  affects  a  cell  within  an  organ  in 
accordance  with  a  simple  Poisson  process  with  rate  parameter  A.  That  is,  the 
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chance  of  an  uninitiated  cell  being  initiated  in  time  interval  (f,  t+h)  is 
approximately  Xh.  If  a  cell  is  initiated  during  exposure  time  it  clones  itself  into 
other  cells  at  rate  p;  the  original  cells  and  its  clones  die  randomly  at  rate  5.  All 
cells  in  the  organ  perform  thus  independently,  according  to  the  model. 
Depending  upon  the  values  of  and  5  (birth  and  death  rates  respectively)  a 
colony  of  initiated  cells  (pre-cancerous,  presumably)  either  tends  to  grow 
exponentially,  or  to  die  off  to  zero  (also  exponentially  fast).  The  fates  of 
colonies  characterized  by  the  same  values  of  birth  rate  and  death  rate  may 
actually  be  entirely  different,  as  befits  experience  with  variability  characteristic 
of  the  real  biological  world.  This  behavior  is  roughly  analogous  to  that  of  the 
flipping  of  the  same  coin:  on  one  occasion  10  flips  may  well  result  in  an  excess 
of  5  Heads  (7  Heads  and  3  Tails),  analogous  to  more  births  (Heads)  than  deaths 
(Tails);  on  another  sequence  of  10  flips  with  the  same  coin  the  result  may  be 
exactly  reversed  (7  Tails,  3  Heads).  Processes  analogous  to  coin  flipping  or  dice 
rolling  can  describe  much,  but  possibly  not  all  interesting  biological  variability 
pertinent  to  risk  analysis.  Other  options  are  suggested  later. 

The  values  of  P  and  <5  describe  clone  colony  properties  in  a  precise 
probabilistic  manner  if  the  model  is  correct.  It  is  only  certainly  approximate, 
but  may  still  provide  a  useful  tool  for  quantifying  risk  of  tumor  formation. 
The  second  step  in  the  malignant  cell  development  process  is  postulated  to  be 
promotion.  A  model  for  this  is  that  at  rate  ji,  i.e.  with  probability  ph  in  time 
it,  t+h)  a  promotion  event  occurs  that  affects  one  of  the  clone  colony  members 
in  proportion  to  the  current  size  of  the  colony;  such  events  are  assumed  to 
occur  in  accordance  with  a  Poisson  process  with  rate  proportional  to 
instantaneous  clone  population  size.  At  the  instant  that  the  first  such 
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promotion  event  occurs  the  clone  colony  (if  one  exists,  i.e.  has  been  initiated) 
will  be  said  to  have  developed  a  tumor,  at  least  in  informal  layman's  terms. 
Note  that  all  original  cells  in  an  organ  are  assumed  to  be  independently 
exposed  to  initiation  and,  thereafter,  to  promotion.  Therefore  all  organ  cells 
and  subsequent  clones,  if  any,  must  survive  from  initiation  to  the  end  of  the 
observation  period  without  being  promoted  in  order  for  the  organ  to  survive 
throughout. 

The  probabilistic  mechainism  described  has  been  used  to  obtain  a  formula 
for  the  survival  probability  for  an  organ  for  any  observation  time  t.  See 
Appendix  A  for  the  formula  and  its  derivation.  Similar  formulas  have’  been 
derived  also  by  Moolgavkar  and  others.  Our  formula  provides  the  basis  for 
statistically  estimating  from  pathology  data,  (combinations  of)  the  parameters: 
X,  the  initiation  rate;  }x,  the  promotion  rate;  and  13  and  5,  the  clonal  birth  and 
death  rates.  Such  estimates  can,  in  turn,  be  used  to  estimate  the  probability  of 
cell,  and  organ,  survived  for  any  time  period.  We  also  supply  an  alternative 
model,  based  on  diffusion  theory  (use  of  mathematical  Brownian  motion)  that 
adds  certain  intuitive  features.  Appendix  B  contains  a  discussion  of  maximum 
likelihood  estimation  from  data  so  as  to  specify  parameters  of  a  preliminary 
model.  Appendix  C  uses  the  preliminary  statistical  model  of  Appendix  B  to 
describe  a  particular  data  set.  Further  work  is  required  to  obtain  additional 
statistical  models  and  procedures  to  analyze  other  experimental  data. 

(d)  Extensions  to  the  Model:  Extra-variation  of  Parameters 

The  above  model,  and  the  conseqxiences  thereof  in  the  form  of  a  survival 
probability  function,  are  appealing  since  they  have  some  plausible  biologic  J 
basis  and  represent  a  form  of  variability  from  organ  to  organ  outcomes 
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(malignancy  evident,  or  not)  that  is  anticipated.  On  the  other  hand,  variability 
between  organs  in  different  subjects  (e.g.  fish)  is  not  well  represented. 
Different,  but  superficially  identical,  biological  entities,  be  they  fish,  rats,  or 
humans,  can  be  expected  to  have  some  differences;  specifically  these  may  be 
between  the  effective  parameters  X,  p.  P  and  5.  If  the  above  are  estimated  from 
data  without  recognizing  the  possibility  of  extra-variation,  biased  results  will  be 
obtained.  See  Harris  (1990)  for  biological  explanations  of  inter-organ  (subject) 
variability. 

There  are  two  simple  and  preliminary  ways  of  dealing  with  the  above 
problem.  One  is  by  attempting  to  "explain"  parameter  variation  by 
representing  it  as  a  function  of  some  causal  variable,  such  as  the  age,  sex, 
weight,  etc.,  of  the  host  subject.  The  technique  could  be  some  form  of 
regression  analysis;  methods  of  McCullagh  and  Nelder  (1983)  suggest 
themselves.  A  description  of  a  preliminary  computational  procedure  to 
estimate  model  parameters  is  described  in  Appei.dix  B.  This  procedure  is  used 
to  estimate  model  parameters  for  a  particular  data  set  in  Appendix  C.  A  second 
approach  is  to  assume  that  the  variability  between  individual  host  organs  can 
be  represented  by  treating  some  or  all  of  the  parameters  as  random  variables 
with  their  own  distributions.  A  typical  survival  function  is  then  obtained  by 
mixing:  the  parametric  survival  function  of  Appendix  A  is  "simply" 
randomized  according  to  the  (joint)  distribution  of  the  parameters.  In  principle 
it  is  desirable  to  recognize  both  sources  of  variability  between  individuals, 
adjusting  for  known  sources  of  variation  by  a  regression  technique  where 
possible,  but  recognizing  the  "unexplainable"  variation  by  use  of  a  mixing 


technique.  The  latter  has  been  carried  out  to  a  limited  extent,  see  Gaver  and 
Jacobs  (1992). 

(e)  Recommendation 

Further  work,  both  theoretical  and  applied,  is  required  to  put  the  above 
ideas  for  characterizing  and  explaining  extra  variability  into  practice.  It  is 
proposed  that  this  work,  plus  effort  to  characterize  the  errors  or  uncertainties  of 
survival  probabilities,  be  carried  out  in  a  follow-on  to  the  current  project. 

3,  CONTACT  WITH  RELATED  RESEARCHERS  AND  INSTITUTIONS 

We  have  established  information-transfer  and  possible  collaborative 
relations  with  several  other  establishments  having  related  objectives.  These 
are 

(a)  The  National  Institute  of  Environmental  Health  Sciences  (NIEHS), 
Research  Triangle,  NOCAR 

Specific  contacts  have  been  made  with  David  Hoel,  Joseph  Hasen:an,  and 
Chris  Portier.  They  have  forwarded  papers.  Hoel  has  invited  me  to  spend 
some  time  there,  which  I  plan  to  do  in  the  spring  of  FY92. 

(b)  RAND  Corp.,  Santa  Monica,  CA, 

A  group  interested  in  environmental  epidemiology,  /'.mong  them  are 
Naihua  Duan  and  Sandy  Ceschwind. 

(c)  EP A,  Washington,  DC 

A  Group  interested  in  risk  analysis.  Dr.  Herman  Gibb. 

(d)  Naval  Medical  Research  Institute,  Toxicology  Detachment,  Wright- 
Patterson  AFB,  Dayton,  OH. 

Particularly  Dr.  Robert  Carpenter. 


6 


We  also  explore  contacts  with  Dr.  Alice  Whittemore,  Professor  of 
Epidemiology  and  Biostatistics  at  Stanford  Uriv.,  and  with  Prof.  Nicholas 
Jewell  of  UC  Berkeley,  and  elsewhere  in  academia. 
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APPENDIX  A.  TWO-STAGE  MODELS 


In  this  appendix  we  present  two  models  for  the  distribution  of  time  until  a 
normal  cell  becomes  promoted.  Tne  first  model  is  a  birth-and-death  model. 
The  second  model  is  a  diffusion  model. 

1.  A  BIRTH-AND-DEATH  MODEL  FOR  THE  TIME  UNTIL  A  NORMAL 
CELL  BECOMES  MALIGNANT  (IS  PROMOTED) 

We  first  develop  an  expression  for  the  distribution  of  dme  until  an 
initiated  cell  or  one  of  its  descendants  becomes  malignant. 

Assume  that  there  is  one  initiated  cell  at  time  0.  Such  cells  divide  at  an 
exponential  rate  P,  and  die  at  an  exponential  rate  5.  Any  initiated  cell  turns 
malignant  at  an  exponential  rate  fi;  i.e.  u  is  the  promotion  rate. 

(a)  Time  to  Promotion  of  an  Initiated  Cell 

Let  T  be  the  random  time  at  which  some  initiated  cell  or  its  descendent 
turns  malignant;  note  that  T  may  actually  be  infinite  if  the  population  of 
initiated  cell  and  its  descendants  dies  out.  Put 

z(t)  =  P{T>t). 

Then  simple  probability  arguments  give  the  equation 

=  g-{p+S+n)t  ^ + 

Q  +  S  +  ul  J  ■' 

^  ^  0 
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Differentiating  with  respect  to  t  and  simplifying  gives 


-z(f)  =  +  S  +  z{t)  - 


=  ^{tf-{p^5+ti)zit)  +  5. 

Hence  2(0  satisfies  a  Riccati  equation  with  initial  condition 

2(0)  3  1 


The  solution  to  (A.2)  with  initial  condition  (A.3)  is 


2(0  = 


where  are  the  solutions  to  the  quadratic  equation 

if  S  3 

x+-  =  0; 

\  P  P  j  P 


if  3  a]  f  3  a]  5  ^ 

p,2=-  1+-  +  -  ±  1  +  -+-  -4-  . 

2[  P  P  [  P  P  P 


P2  5  1  and 


Since  J  - 


1  i 
1  +  0  +  ^ 


4^j  ^  1  +  ^  -r-  ^  ,  both  pj  and  P2  are  positive.  Further 


f  5  uY  52 
Pl-P2=  +  ^  +  ^  -4-^  >0- 
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Hence, 


lim  2(0  =  P2- 

If  the  death  rate  5=0,  then  p-,  =  0  and  lim  P(T  >  t)  =  0;  if  5  »  0,  then  there  is  no 

*■  i 

death  of  initiated  cells  and  thus  an  initiated  cell  will  transition  to  a  malignant 

ceil  in  a  finite  time  with  probability  1.  If  5>  0,  then  the  initiating  cells  can  die, 
thus  preventing  a  transition  to  malignancy  and  hence  ^lim  PIT  >  f)  =  P2  >  0. 


(b)  A  Model  for  the  Time  until  a  Normal  Cell  becomes  Malignant  (is 
promoted) 

Assume  that  each  normal  cell  is  initiated  at  an  exponential  rate  Ag.  Let  N 

% 

be  the  total  number  of  normal  cells  in  an  organ.  Let  S  denote  the  first  time  a 
normal  cell  transitions  to  a  malignant  ced. 


P{SSf}  = 


(A.7) 


where  z  is  a  given  in  (A.4).  Assume  Ag  is  small  and  put  A  =  XqN,  a  constant. 
Then 


P(S  >  t}  »  exp 


Nlnl 


(A.8) 


‘  exp-j  -  At  +  X  j  z{s)ds  I 

0  J 


(A.9) 


exp 


X(p^-l)t-X~\n\ 


P1-P2 


(A.IO) 
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2.  A  DIFFUSION  APPROXIMATION  TO  THE  TIME  TO  PROMOTION 

An  alternative  approac±i  to  constructing  a  model  for  time  to  promotion,  i.e. 
appearance  of  a  tumor,  is  via  diffusion  approximation;  see  Karlin  and  Taylor 
(1981),  Chap.  15.  The  approach  can  be  mathematically  justified  as  an  asymptotic 
limit  of  a  sequence  of  discrete-state  processes  by  use  of  weak  convergence 
theory;  see  Kopp-Schneider,  Portier  and  Rippman  (1991)  whc  have  used  this 
methodology  also  and  discussed  the  weak  convergence  issues.  Here  our 
motivation  is  to  model  biological  phenomena,  i.e.  variability,  ir.  a  convenient 
and  flexible  manner,  so  questions  of  mathemafral  rigor  can  be  temporarily 
suppressed. 

Consider  'n  organ  that  at  time  0  contains  x  cells  that  have  just  been 
initiated.  Thiiik  of  x  as  a  real  fixed  positive  number,  not  insisting  that  it  be  an 
integer.  Now  suppose  promotion  is  a  Poisson  process  of  rate  n,  applying 
independently  to  each  of  the  x  cells.  This  impl*  i  that  no  promotion  occurs  in 
time  (f,  f+/j)  with  probability  »  l-fih  *  o{h)  is  h~*0.  The  diffusion 
approximation  to  the  growth/death  of  the  entire  initiated  cell  population,  of 
size  /(O)  ■  X  initially,  changes  to 

lih)  a  X  + 

at  time  h  Ox  small  after  initiation),  with 

tx  Normally  distributed  with  mean  v{x)h  oih),  and  variance  (7*(x)/«+o(/i) 
Example.  A  diffusion  approximation  to  the  previous  birth-death  model  would 
naturally  equate  the  mean  and  variance  of  change  in  a  small  time  interval  of 
length  h  to  v(x)h  and  cr^(x)/i  respectively;  these  turn  out  to  be 

u(x)/t  =  x(/J-o)/i  (A.U) 
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and 


o^(x)h  =>  x(p+5)k  (A,12) 

and  the  birth-death  model's  increment  distribution  (net  births  over  deaths)  by  a 
normal  variable  with  the  matched  moments. 

The  parameters  tXi)  and  (P'i.x)  are  known  as,  respectively, 
iXx):  Infinitesimal  mean,  or  drift; 

o^Cx);  infinitesimal  variance,  or  diffusion  coefficient. 

The  u(x)  parameter  describes  the  tendency  of  the  process  to  move 
deterministically  in  a  given  direction  during  a  time  of  length  h;  the  ff^ix) 
parameter  provides  a  scale  ^or  the  random  variability  around  the  mean. 

Put  zU;  x)  for  the  survival  probability,  until  tu  nor  apf?earance,  given  that  x 
cells  are  initially  initiated: 

zit;  x)  -  P{T,  >  t)  -  P(T  >  1 1  /(O)  -  x). 

The  assumptions  made  initially  then  suggest  that  z  satisfies  the  following 
equation  if  «  1: 

z(t;x)a  \e~^- - - - - - zU-h.x-^u).  (A.13) 

'/Zx^jcr{x)h 


Now  as  h  is  small  then  y  will  also  assume  small  values  with  appreciable 
probability  so  replace  zU-h,  x+y)  by  its  Taylor  series  expansion 


,,  ,  \  \  kdz  dz  v~ 


Then  the  integral  can  be  performed  and 
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13 


(A.18) 


-  ^2)  -  ‘52(1  -  ^l)exp|^Vu^  +  2M(T^}} 
(1  -  <52)  -  (1  -  4i)exp|^V^7+2^?  jlj 


where 


,  V  .  1 


(A.19) 


Thus,  if  at  time  0  there  is  one  initiated  cell,  then  the  probability  that  this 
cell  or  one  of  its  descendants  is  not  promoted  by  time  t  is 

P{T  >  t)  =  exp[(jp(t)}  (A.201 

where  <p(t)  is  given  by  (A.18).  This  expression  is  to  be  compared  to  the  birth* 
and-death  model  expression  (A  4).  Note  that  as  t-*—,  PIT  >  t)  exp{(®2),  which 
is  always  positive.  A  further  argument  similar  to  that  in  (A.7)  is  needed  to 
obtain  the  distribution  of  time  until  a  normal  cell  Is  initiated  and  promoted. 
Further  discussion  of  the  relationship  between  these  models  will  be  given 
elsewhere. 
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APPENDIX  B.  MODEL  HTTING  METHODS  AND  QUANTIFYING  BIOASSAY 

SURVIVAL 

L  PRELIMINARY  STATISTICAL  MODELS  AND  METh  UDS  FOR 
ANALYZING  BIOASSAY  DATA 

Suppose  N  organisms  (for  example  fish)  are  used  in  an  experiment. 
Groups  of  these  organisms  may  be  exposed  to  different  treatments.  Let  T-  be 
the  random  time  until  organism  i  develops  a  particular  symptom,  e.g.,  cystic 
degeneration.  Let  =  (X^,  ...»  X,-^  be  covariates  which  (possibly)  influence 

Tj;  the  Xj  could  be  levels  of  substances  having  possible  toxic  effects  to  which 
the  organisms  are  exposed.  Let  G(t;  xp  »  P[Tj  5 1 !  X^  »  xj.  We  will  assume  that 
the  organisms  develop  symptoms  independently  of  each  other.  In  this  initial 
model,  the  symptom  is  either  present  or  not. 

Suppose  that  organisms  are  sacrificed  at  time  with  <  fj  <  •••  <  • 

We  will  label  the  organisms  so  that  organisms  1  through  are  sacrificed  at 
time  organisms  «j+l,  ...,  rtj+n  j  are  sacrificed  at  time  t2;etc.  Let  S;  =»  1  if 
organism  i  exhibits  the  symptom  when  it  is  examined.  Under  the  assumption 
of  independence,  the  likelihood  function  is 

Kny  __ 

*>.1  1-1 

where  ng  =  0  and  G(t;  x)  =  \-C{t;  x).  The  likelihood  functions  form  the  basis 
for  estimation  of  param.eters  in  the  distributions  that  model  survival  times,  i.e. 
G. 
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Example  (Simple  Binomial  Model).  If  there  are  no  covariates,  then  (B.l) 
becomes 


L=n(A)c('k/*G(<»)''‘'^  (B.2) 

k^\ 

where  is  the  number  of  the  organisms  exhibiting  the  symptom. 

A  procedure  to  estimate  the  parameters  of  the  distribution  G  for  the  simple 
binomial  model  is  as  follows. 


2.  MAXIMUM  LIKELIHOOD  ESTIMATION  IN  THE  SIMPLE  BINOMIAL 
MODEL 

(a)  l  ikelihood  and  Parameter  Estimation  Formulas 

Assume  the  distribution  of  the  time  to  appearance  of  a  symptom,  G,  is  a 
function  of  the  parameters  aj ,  /  =  1, In  this  section  we  discuss  maximum 
likelihood  estimation  of  for  the  simple  binomial  model.  Presumably  the 
subjects  examined  at  time  t,^k  =  1,2, K  have  all  been  subjected  to  a  common 
dosage  of  a  potential  toxin.  The  purpose  of  the  present  analysis  is  to  predict 
survival  probabilities  as  they  depend  on  such  dosage.  The  log-likelihood 
function  for  the  simple  binomial  model  is 

1=  )  +  /klnG(f)t;a)  +  (nfc-/i.)lnG(ft;a)  (B.3) 

1:=1 


where  a  =  ctj).  Differentiating,  we  obtain 


fk 


ktiG(fk;ct)  aay 


G(t^;a) 


C((*;cr) 
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J:=l 

K 

I 

ib-1 


Since  E[f  J  =  ny^G{t,,;  a) 
3^ 


G(f]fc;a)G(fjt;a) 

■G{tk;a). 


—G(t,,;a) 

; 


fk-nkG{tk;a) 

3 

G{tfc;a)G(tit;a)_ 

daj 

(B.4) 


daida„ 


-I-t-' - = - 


]k=l 


G(ijt;a)G(fit;a) 


(B.5) 


Thus  a  Newton  procedure  for  finding  the  maximum  likelihood  estimates  of 
lap  j  =  1, ...,/}  would  iteratively  solve  the  system  of  linear  equations 


’’  ^  ^ 


m-1  J 


(a, -a?) 


(B.6) 


where  a  = 


f  0  0^ 

a^,...,aj 


Such  iterative  procedures  can  be  programmed  for  a 
digital  compiiter,  and  the  resulting  parameter  values  can  be  used  to  compute 
predictions  for  survival  probabilities,  or  risk,  as  the  latter  depend  upon  the 
parameters  of  such  models  as  described  in  Appendix  A.  An  illustration  is 
given  in  Appendix  C. 


APPENDIX  C  A  PRELIMINARY  ANALYSIS  OF  DATA 


We  present  and  disooss  a  specific  data  analysis  and  uncertainty  assessment 
utilizing  the  methodology  described  above.  Note  that  this  analysis  is  less 
comprehensive  than  that  potentially  possible. 

The  data  analyzed  are  taken  from  the  U.S.  Army  Biomedical  Research  and 
Development  Laboratory  Study  N;  Utilization  of  Fish  to  Evaluate  Carcinogenic 
Potential  of  Tricholoethylene  Contamined  Groundwater — Pathology  Report, 
Vol.  1  of  3.  The  report  presents  results  of  a  "histopathologic  examination  of 
tissues  from  fish  of  the  spedes  Oryzias  latipes  (Japanese  medaka)  which  was 
performed  to  evaluate  the  cardnogenic  potential  of  trichloroethylene  (TCE)  in 
groundwater.  Fish  exposed  to  groundwater  with  various  concentrations  of 
TCE  were  either  pre-treated  or  not  pre-treated  with  10  mg/1  diethylnitrosamine 
(DEN)  in  water  for  48  hours  at  17  days  post  hatch.  Exposures  were  begun  in  the 
biomonitoring  trailer  on  the  sixth  day  after  treatment  with  DEN."  After  three 
months  of  exposure  25  fish  in  each  treatment  group  were  sacrificed,  (interim 
sacrifice).  The  remaining  fish  in  each  of  the  treatment  regimes  were  divided 
into  two  groups  "one  group  of  fish,  designated  the  chronic  fish,  continued  to  be 
exposed  to  various  concentrations  of  TCE  in  the  water  for  an  additional  three 
months."  Tne  chronic  fish  were  sacrificed  at  6  months  into  the  study.  It  is  the 
data  on  these  fish  that  we  will  consider. 

The  data  we  consider  are  the  occurrence  or  nonoccurrence  of  cystic 
degeneration  (CD)  in  the  chronic  fish  in  each  treatment  group.  The  data  are 
presented  in  Table  1. 
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TABLE  1. 


#  fish  with  symptom/#  fish 
killed 

Sacrifice  Time 


%TCE 

3  months 

6  months 

0 

6/25 

4/15 

25 

4/25 

5/13 

50 

2/25 

4/14 

100 

3/25 

3/14 

0 

11/25 

6/12 

25 

4/25 

8/13 

50 

6/25 

5/12 

100 

7/25  i 

3/8 

The  statistical  model  used  is  the  simple  binomial  model  of  Appendix  B 
with  exponential  distribution 

G(t)  =»  l-€xp{-e>?). 

The  parameter  -y  is  the  natural  logarithm  of  mean  time  to  exhibit  symptoms. 
The  Newton  procedure  described  in  Appendbc  B  was  used  to  estimate  the 
parameter  y  numerically  for  the  data  in  each  group.  Each  group  contains  two 
data  points,  one  for  three  months  and  one  for  six  months.  The  maximum 
likelihood  estimates  appear  in  Table  2. 


TABLE  2.  MAXIMUM  LIKELIHOOD  ESTIMATES  OF  MEAN  TIME  TO 

EXHIBIT  CD 


Group 


2 

X 

4 

X 

6 

X 

8 

X 

%TCE 

-y 

=  Estimated  Mean 
Time  to  CD 

0 

2.66 

14,3 

25 

2.68 

14.5 

50 

!  3.32 

23.8 

100 

3.28 

24.4 

0 

1.85 

6.4 

25 

2.31 

10.0 

50 

2.40 

11.0 

100 

2.32 

10.2 

The  last  column  of  Table  2  shows  that  the  estimates  of  the  mean  time  to  CD  for 
those  groups  in  which  the  fish  were  pretreated  with  DEN  are  always 
considerably  smaller  than  those  for  fish  that  were  not  pretreated.  The 
estimates  are  based  on  use  of  a  very  simple  model,  the  exponential.  Note, 
however,  that  the  two-stage  clonal  expansion  model  of  survival  probability 
described  in  Appendix  A  has  approximate  exponential  behavior  in  the  right 
tail,  so  choice  of  the  exponential  for  data  analysis  is  not  inconsistent  with  that 
theory. 

There  has  been  no  opportunity  to  assess  possible  between-fish  variability  in 
susceptibility  to  incur  different  possible  toxin  effects.  This  is  for  the  future. 

VARIABILITY  OR  UNCERTAINTY  ASSESSMENT 

In  order  to  assess  uncertainty  in  the  estimate  of  7  and  the  mean  time  to  CD 
appearance  the  technique  of  bootstrapping  was  used,  cf.  Efron  and  Tibshirani 
[19S6].  This  is  a  re-sampling  methodology  that  is  of  wide  and  useful 
application.  For  each  group  treatment,  experimenlal  data  was  simulated  using 
the  estimated  model  using  the  y-estimate  in  Table  2;  that  is,  for  each  bootstrap 
replication  for  each  group  of  fish,  two  independent  binomial  random  variables 
were  simulated;  one  with  25  trials  and  probability  of  success  l-exp{-3e"V),  and 
the  other  having  probability  of  success  1-exp {-6e“V)  with  the  number  of  trials 
equal  to  the  number  of  fish  in  the  chronic  population  sacrificed  at  six  months. 
Any  fish  that  have  presumably  died  (or  vanished)  were  simply  ignored;  note 
that  these  could  have  been  removed  because  they  very  early  reached  a  dose- 
related  fatal  endpoint,  and  hence  the  treatment  of  these  is  a  potential  source  of 
bias.  These  simulated  data  were  then  used  to  obtain  one  bootstrap  estimate  of 
7.  For  each  treatment  group  500  independent  bootstrap  estimates  of  7 were 
obtained. 
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Figure  1  presents  boxplots  of  the  bootstrap  estimates  of  -7  for  each 
treatment  group;  recall  that  -7  corresponds  to  the  logarithm  of  the  mean  time 
to  symptom  development. 

The  following  is  taken  from  the  documentation  of  GRAFSTAT,  a 
developmental  product  of  IBM  which  the  Naval  Postgraduate  School  is  using 
under  a  test  agreement  with  IBM.  'The  box  portion  of  the  plot  extends  from 
the  lov/er  quartile  of  the  sample  to  the  upper  quartile.  (The  lower  quartile  is 
the  point  for  which  one  quarter  of  the  sample  lies  below  and  three  quarters 
above.  The  upper  quartile  is  analogous.)  The  line  across  the  center  of  the  box 
marks  the  median.  The  circle  in  the  box  represents  the  mean. 

The  distance  from  the  lower  to  the  upper  quartile  is  called  the  interquartile 
distance,  and  it  will  be  represented  by  Q.  The  points  at  the  ends  of  the  two  lines 
(called  whiskers)  are  the  smallest  and  largest  points,  respectively,  within  1.5Q 
of  the  quantiles.  The  points  beyond  the  whiskers  arc  outlying  values." 

The  boxplots  for  groups  2,  4,  6  and  8  (those  groups  pretreated  with  DEN) 
generally  lie  below  those  for  groups  1,  3,  5,  and  7  respectively  (those  groups  not 
pretreated  with  DEN).  Thus  treatment  with  DEN  tends  to  shorten  the  mean 
time  to  symptom  development,  even  accounting  for  sampling  errors  in  the 
estimates.  The  bootstrapping  results  tend  to  strengthen  inference  concerning 
the  effect  of  DEN  treatment,  and  they  provide  perspective  on  the  sensitivity  of 
the  experimental  procedures,  e.g.,  how  well  dosage  effects  can  be  distinguished. 
Comparison  of  the  width  of  the  boxes  of  groups  2,  4,  6,  and  8  with  those  of  1,  3, 
5,  and  7  suggests  that  pretreatment  with  DEN  tends  to  decrease  the  amount  of 
variability  of  the  y-estimates  (as  well  as  reducing  those  estimates).  The 
variability  of  the  estimate  of  yprecludes  any  strong  conclusions  concerning  the 
effect  of  the  different  dosages  of  TCE  on  corresponding  mean  times  to  symptom 
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development.  The  design  of  the  experiment  needs  to  be  modified  to  increase 
the  experiment's  sensitivity  if  higher  sensitivity  is  required. 
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